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We provide an exact analytical solution of the collapse dynamics of self-gravitating Brownian 
particles and bacterial populations at zero temperature. These systems are described by the 
Smoluchowski-Poisson system or Keller-Segel model in which the diffusion term is neglected. As 
a result, the dynamics is purely deterministic. A cold system undergoes a gravitational collapse 
leading to a finite time singularity: the central density increases and becomes infinite in a finite 
time t co u. The evolution continues in the post collapse regime. A Dirac peak emerges, grows and 
finally captures all the mass in a finite time t en d, while the central density excluding the Dirac peak 
progressively decreases. Close to the collapse time, the pre and post collapse evolutions are self- 
similar. Interestingly, if one starts from a parabolic density profile, one obtains an exact analytical 
solution that describes the whole collapse dynamics, from the initial time to the end, and accounts 
for non self-similar corrections that were neglected in previous works. Our results have possible 
application in different areas including astrophysics, chemotaxis, colloids and nanoscience. 



I. INTRODUCTION 

In a series of papers (seeQ for a short review), we have 
studied a class of mean field drift-diffusion equations of 
the form 



dt 



- (Vp + pV$) 



'at 



A$ - fc 2 $ - S d G P , 



(1) 



(2) 



where the equation of state p = p{p) can take different 
shapes. These equations, derived in a statistical mechan- 
ics context in 0, can be viewed as nonlinear mean 
field Fokker-Planck (NFP) equations [J] describing a sys- 
tem of Langevin particles in interaction f6lj . The first 
equation (JXJ) can be interpreted as a generalized Smolu- 
chowski equation in which the evolution of the density 
p(r, t) is due to a competition between a (nonlinear) diffu- 
sion and a drift. Contrary to the ordinary Smoluchowski 
equation the potential $(r, t) is not fixed but de- 
termined by the density itself according to the reaction- 
diffusion equation (|2|). Equations (d])-© are based on 
a mean field approximation which is known to be exact 
for long-range interactions when the number of particles 
N — > +00 @ ■ On the other hand, the linear diffusion 
V • (DV/o) in the ordinary Smoluchowski equation is re- 
placed by a more general term V- (£ _1 Vp) which can lead 
to anomalous diffusion. The function p(r, t) = p[p(j, t)] 
can be interpreted as a barotropic pressure It can 
take into account microscopic constraints such as short- 
range interactions, close packing effects, steric hindrance, 
non-extensivity, exclusion or inclusion principles... Ac- 
cordingly, the drift-diffusion equations (HJ-© are asso- 
ciated with generalized forms of free energy that can be 
non-Boltzmannian. This leads to a notion of effective 
generalized thermodynamics (EGT) 0, 0] ■ 



The drift-diffusion equations ([I])-© can be derived 
heuristically from generalized stochastic processes in 
physical space of the form 




R(t), 



(3) 



where $(r, t) is a mean field potential determined by Eq. 
© and R(t) is a Gaussian white noise satisfying (R(t)) = 
and (Ri(t)Rj(t')) = %<5(t - t') where i = 1, ...,d label 
the coordinates of space. The Fokker-Planck equation as- 
sociated with the generalized Langevin equation ([3]) is the 
generalized mean field Smoluchowski equation (JlJ. Note 
that the strength of the noise can depend on the local 
density p(r, t) of particles which may lead to anomalous 
diffusion [8j. The generalized mean field Smoluchowski 
equation (|T|) can also be derived from the master equation 
by assuming that the probabilities of transition explicitly 
depend on the occupation numbers of the initial and ar- 
rival states (see @ and Sec. 2.11 of (U). This kinetical in- 
teraction principle (KIP) takes into account microscopic 
constraints such as exclusion or inclusion principles |9|. 

The model (HJ-© assumes an overdamped evolution in 
which the velocity of the particles is directly proportional 
to the force — V$ as in Eq. ((3]), but more general mod- 
els taking into account inertial effects can be introduced 
@. They are based on generalized stochastic processes 
in phase space of the form 



dr 

dt 



(4) 



dv 



l2Df 



C(f) 



f 



R(t), (5) 



where the noise depends on the distribution function 
/(r,v,t) of the particles. Here, £ is a friction coeffi- 
cient and C(f) is a convex function (C" > 0). This 
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function determines a generalized entropic functional 
S = — J C(f) drdv The Fokker-Planck equation cor- 
responding to the generalized Langcvin equation (0 is 
the generalized mean field Kramers equation 



to the generalized Darcy law 



df 
dt 



— - V$ ■ — 
dr <9v 



d_ 



TC"(f)f^ + /v 



e^J- = A$ - fc 2 $ - S d Gp, 
at 



(6) 
(7) 



where we have defined the generalized thermodynamical 
temperature T through an Einstein-like formula D = £T 
[H . These equations govern the evolution of the distribu- 
tion function /(r, v, i) in phase space. The generalized 
Smoluchowski equation ([!} can be derived from the gen- 
eralized Kramers equation ([5]) in a strong friction limit 
£ — > +00 by using a Chapman-Enskog expansion 
or a method of moments [3. In the strong friction limit, 
there exists a precise link [3[ between the barotropic pres- 
sure p(p) occuring in the generalized Smoluchowski equa- 
tion |T]) and the function C(f) occurring in the gener- 
alized Kramers equation (J6j) . For example, the Boltz- 
mann entropy Sb = ~ / m ^ n m ^ r ^ v (leading to normal 
diffusion) |5| determines an isothermal equation of state 
p(r,t) = p(r,t)kBT/m, where T is the temperature and 
m the individual mass of the particles. The Tsallis en- 
tropy S q = —jzj J (f 9 — f) drdv (leading to anomalous 
algebraic diffusion) [ll[ determines a polytropic equa- 
tion of state p(r,t) = Kp(r,t) 7 , where K is the poly- 
tropic constant and 7 the polytropic index. They can 
be related to T and q as explained in 0). We can also 
consider Fermi-Dirac and Bose-Einstein entropies (taking 
into account exclusion or inclusion constraints)determin- 
ing fcrmionic and bosonic equations of state |3[. For ex- 
ample, the equation of state p{p) = — Too ln(l — p/ao) 
that takes into account an exclusion constraint in physi- 
cal space has been studied in [l2j . 

The generalized Smoluchowski equation ((Tl can also 
be derived from the damped Euler equations [2, El : 

^+V-(pu) = 0, (8) 

-^ + (u-V)u = — Vp-V$-£u. (9) 
at p 

The first equation is the equation of continuity (taking 
into account the local conservation of mass) and the sec- 
ond equation is the momentum equation. We have as- 
sumed that the pressure is isotropic and barotropic, i.e. 
it is a function of the density p = p{p). In the ideal case, 
it is given by the perfect gas law p(r, t) = p(r, t)kBT/m. 
More generally, we allow p{p) to be nonlinear so as to take 
into account microscopic constraints such as short-range 
interactions and other non ideal effects. In the strong 
friction limit £ — > +00, we can formally ignore the iner- 
tia of the particles in the damped Euler equation leading 



1. 



-Vp- V$. 



(10) 



Substituting this relation in the continuity equation ([8]). 
we obtain the generalized Smoluchowski equation (fTj). 
We note, however, that the damped Euler equations ((HJ- 
© cannot be directly derived from the stochastic equa- 
tions (01)- © or from the Kramers equation (J6j> as they 
rely on a local thermodynamic equilibrium (LTE) assump- 
tion that is not rigorously justified. By contrast, the 
generalized Smoluchowski equation (fTj) can be rigorously 
derived from the generalized Kramers equation ( l 6 | ) in a 
strong friction limit £ — > +00 as shown in Refs [3|, [H| . 

We finally note that the generalized mean field Smolu- 
chowski equation ([T l)-([^1) has recently been derived from 
a kinetic theory [13( that combines the technics used in 
the theory of systems with long-range interactions and in 
the theory of liquids (see Appendix |C|) . In this approach, 
the long-range interactions (leading to the drift) are mod- 
eled by using the mean field theory and the short-range 
interactions (leading to the pressure term) are modeled 
by using the dynamic density functional theory (DDFT) . 
This new derivation (l3j | provides an alternative to the ki- 
netic theory based on generalized thermodynamics 0, H| . 

The generalized mean field drift-diffusion model ([1])- 
(P| appears in a number of physical situations that we 
briefly review. 

(i) When $ = <fr e xt (r) is an external potential and p = 
pksT '/m, we recover the ordinary Smoluchowski equation 
describing, for example, the sedimentation of colloidal 
particles in a gravitational field Q. When $ e2 . t = and 
p = Kp" 7 , we recover the porous medium equation flil ]. 
The porous medium equation with an external potential 
<fr ext (r) has been considered by Plastino & Plastino [HI 
in connection with Tsallis generalized thermodynamics 
P. 

(ii) When e = k = 0, and for an isothermal equation 
of state p = pkBT/m, we get the Smoluchowski- Poisson 
system. In the repulsive case G < 0, i.e. for the Coulom- 
bian interaction, this corresponds to the Debye-Hiickel 
model of electrolytes [l(|. The same equations, called the 
Nernst-Planck equations [l7j|, are also used to study ion 
transport in biological channels [3, Oil and carrier trans- 
port in semiconductors [2^. In the attractive case G > 0, 
i.e. for the gravitational interaction, the Smoluchowski- 
Poisson system describes a gas of self-gravitating Brown- 
ian particles [l|. In that case, G represents the constant 
of gravity and Sd is the surface of a unit sphere in d 
dimensions (62j . 

(iii) When e > and k > 0, we obtain a generalized 
form of Keller-Scgcl model describing the chemotaxis 
of bacterial populations. In that case, p denotes the den- 
sity of bacteria and c = — $ the concentration of the se- 
creted chemical. Furthermore, the term — fc 2 $ takes into 
account a possible degradation of the chemical, leading 
to a shielding of the interaction on a scale jl2| . 

(iv) Recently, it has been shown [30[ that Eqs. ©-© 
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with e = 0, k > and d = 2, describe the dynamics of 
colloids at a fluid interface driven by attractive capillary 
interactions. In that context, p is the particle density, <3> 
the interfacial deformation and k~ l the capillary length 
whose finiteness leads to a shielding of the interaction. 
Finally, G is equal to the ratio between the capillary 
monopolc associated to a single particle and the surface 
tension. 

(v) Drift-diffusion equations of the form fTjl-pj). or 
their generalization (see Appendix [X]), have been pro- 
posed to describe phase seg regation in model alloys with 
long-range interactions [3l| and aggregation of finite-size 
particles and directed-assembly processes for nanoscience 

(vi) Nonlinear mean field Fokkcr-Planck equations of 
the form ©-© and ©-(HI) have also been introduced in 
two-dimensional (2D) vortex dynamics and stellar dy- 
namics in order to provide a small-scale parametriza- 
tion of the 2D Euler-Poisson and Vlasov-Poisson systems 
[33l H3] . They describe the rapid formation of quasi sta- 
tionary states (QSS) on the coarse-grained scale repre- 
senting large scale vortices or galaxies [35|. In 2D hy- 
drodynamics, p corresponds to the coarse-grained vor- 
ticity u; and $ to the stream function ip. On the other 
hand, Eq. © reduces to the Poisson equation —Atf,' = uJ. 
In the case of geophysical flows described by the quasi 
geostrophic (QG) equations, the Poisson equation is re- 
placed by a screened Poisson equation — A-0 + -0/i? 2 = q 
where R is the Rossby radius and q the coarse-grained 
potential vorticity. In stellar dynamics, / is the coarse- 
grained distribution function, p is the density of stars and 
<3? the gravitational potential. 

The numerous analogies between these very different 
systems, covering all scales of physics (self-gravitating 
systems, large-scale vortices, biological organisms, elec- 
trolytes, colloids, semiconductors, nanotechnology,...) 
show the importance of studying nonlinear mean field 
Fokker-Planck equations of the form ©-© and ©-© 
[63j . Depending on the equation of state p(p, T) , on the 
value of the temperature T and on the dimension of space 
d, the system can either (i) converge towards an equi- 
librium state, (ii) collapse, or (hi) evaporate. A short 
description of these different regimes in the case of equa- 
tions ([I])-© with e = k = and in the case of isothermal 
and polytropic equations of state has been given in 
It reveals the complexity and the richness of these appar- 
ently simple equations. 

In this paper, we treat in detail a simple case cor- 
responding to e ~ k ~ 0, G > (long-range New- 
tonian interaction) and p = (cold system) [64| . In 
other words, we consider the overdamped dynamics of 
particles in Newtonian interaction in the zero temper- 
ature limit T = 0. This system undergoes a gravita- 
tional collapse resulting ultimately in the formation of a 
Dirac peak containing all the mass. As discussed above, 
this result can have application in other areas such as 
chemotaxis, colloid dynamics and nanotechnology. In- 
terestingly, if we start from a parabolic density profile, 



we can obtain an exact analytical solution of the equa- 
tions for all times. This solution describes both pre and 
post collapse regimes. The pre-collapse regime leads to 
the formation of a finite time singularity: the central 
density increases like (t co u — t)^ 1 and becomes infinite 
in a finite time t co u resulting in a singular density pro- 
file proportional to r ~ 2d /( d + 2 ) , The evolution contin- 
ues in the post-collapse regime with the formation of 
a Dirac peak at r = 0. The mass contained in the 
Dirac peak grows like (i co z ; A) (d+2)/2 (£ - t coU ) d / 2 while 
the central density excluding the Dirac peak decreases 
like (t — t con) {t en d/ 't — 1). Finally, in a finite time 
tend = ((d + 2)/2)t co u, all the particles are ultimately 
absorbed by the Dirac peak. In our previous works 
[3ril - [38j . we had only described the self-similar dynam- 
ics of the system at T = close to the collapse time 
t co u and for r — > 0. The present analytical solution 
extends our results to all times < t < t en d and all 
radii < r < r max (t), and describes exact corrections to 
the self-similar solution. It provides a simple illustration 
of the Dirac peak formation in the post-collapse regime. 
These behaviors (pre and post collapse, finite time sin- 
gularity, growth of a Dirac peak) also arise at T ^ 
(with different exponents) but they are more difficult to 
investigate analytically |36l439| . Therefore, the present 
explicit solution can be useful before considering more 
complicated situations. 

II. THE FORMAL SOLUTION OF THE 
PROBLEM 

The Smoluchowski-Poisson (SP) system at T = re- 
duces to the form 

A$ = S d G P . (12) 
Equation (fTTj) is equivalent to a continuity equation 

^+V-(pu)=0, (13) 

with a velocity field 

u(r,t)=-iv*(r,f). (14) 

Introducing a Lagrangian description of motion like in 
[HI , the equation determining the trajectory of a "fluid" 
particle is 

§=u(r,t)=-iv*(r,t). (15) 

It coincides with the equation of motion © when the 
noise is switched off (T = 0). In that case, the dynamics 
is deterministic. If the distribution of particles is initially 
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spherically symmetric, it will remain spherically symmet- 
ric for all times. Using the Gauss theorem, the equation 
of motion (Ti~5l) becomes 



dr _ lGM(r,t) 
dt 



r d-l 



where 



M(r, t) 



(16) 



(17) 



is the mass contained within the sphere of radius r at 
time t. Conversely, 

*'^^¥^ ) - (18) 

Since the equation of motion (|16j) is a first order dif- 
ferential equation, the particles do not cross each other. 
Therefore, the conservation of mass implies 

M(r,t)=M(a,0), (19) 

where r is the position at time t of the particle located at 
r = a at t = 0. Equation (|16[) can therefore be rewritten 

dr _ 1 GM(a, 0) 



and it is easily integrated into 



r d-l 



(20) 



= a d — -GM(a,0)t. 



(21) 



Equations ([191) and ([21]) give the general exact solution of 
the problem [36j . The time at which the particle initially 
located at r = a reaches the origin r = is 



i*(a) = 



(22) 



dGM(a,0) 

This time increases with the distance provided that 
M(a,0)/a d is a decreasing function of a (physical case). 
The opposite situation is treated in Appendix [T~ 

Remark: the gravitational collapse of a cold gas (T = 
0) initially at rest described by the Euler-Poisson system 
has been investigated long ago by Hunter [Z(| , Mestel (4l| 
and Penston [42| in astrophysics. In that case, Eq. (]2T)]) 
is replaced by a second order differential equation which 
can be solved in a parametric form. The collapse of a stel- 
lar system which starts from a configuration in which all 
the stars have zero velocity is initially similar to that of a 
cold gas at T = 0. However, due to fluctuations (finite N 
effects), the velocity dispersion of the stars increases and 
orbit crossing occurs [22[ ■ In that case, the mass interior 
to a given particle changes so that the collapse solution 
becomes much more complicated. The system undergoes 
damped oscillations before finally settl ing on a virialized 
state as a result of violent relaxation [43[ . The simplifi- 
cation obtained by considering self-gravitating Brownian 
particles in the overdamped limit allows us to obtain a 
complete and explicit solution of the collapse dynamics 
for all times since there is no orbit crossing in that case. 
As we shall see, the system creates a Dirac peak at the 
origin instead of exhibiting damped oscillations. 



III. THE COLLAPSE OF A COLD 
HOMOGENEOUS SPHERE 

Let us first consider the collapse of a homogeneous 
sphere of mass M and initial radius R. Its initial density 
is p(a,0) = dM/(S d R d ) for a < R and p(a,Q) = for 
a > R. The corresponding mass profile is 



M(a,0) = ^a d , 



(23) 



for a < R and M(a, 0) = M for a > R. According to 
Eqs. (|2Tj) and (f23|) . the position at time t of a particle 
initially located at r = a is 



r a = 1 - 



dGM 



t a a 



Therefore, the particles reach r = at a time 



£R d 
dGM ' 



(24) 



(25) 



whatever their initial position a. This leads to a complete 
collapse of the system (Dirac peak) at a time t = f*. 
According to Eqs. CED, ([23]) and (J2U), the mass profile 
at time t is 



M{r,t) = 



M 



R d 1 



j_ 

i. 



(26) 



The sphere of particles remains spatially homogeneous 
with a radius decreasing in time like 



R(t) = R 



t 

1 

t * 



l/d 



_ The density increases with time like 



P(t) 



P(0) 



1-T- 



l_ ■ 

t, 



(27) 



(28) 



This result can be directly obtained from Eqs. (fTT]) -(fl"2" j) 
using the fact that the density is spatially homogeneous. 
Indeed, using the Poisson equation (TP2"j) . Eq. (TTT1) be- 
comes 



dp Sd.G 2 



(29) 



which leads to Eq. ([28]) after integration. 

Remark: in the context of the Euler-Poisson system 
at T = 0, the free-fall time of a homogeneous sphere is 
iff = (3tt/32Gpo) 1/2 (in d = 3) [Ij. This is of the order 
of the dynamical time td yn = (37r/16Gpo) 1 ^ 2 [111- F° r 
the Smoluchowski-Poisson system at T = 0, the collapse 
time (|25p of a homogeneous sphere can be written t* = 



i/SdPoG. It is of the order £t 2 



dyn ' 
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IV. THE PRE-COLLAPSE OF A PARABOLIC 
PROFILE 

A. The exact solution 



time singularity in which the central density diverges. 
We note that 



d+2 



(38) 



We now consider an initial mass profile of the form 

M(a,0) = A(a d - Ba d+2 ), (30) 

where A and B are two positive constants. Using Eq. 
(fT8]) . we find that the corresponding density profile is 
parabolic 



dA ( d + 2 , 



(31) 



These expressions are valid for a < R where R is the 
radius at which the density vanishes. The constant B is 
related to the radius R of the initial configuration by 



B 



d 



(d + 2)R 2 



(32) 



On the other hand, the constant A is related to the total 
mass M of the system, and to its radius R, by 



.4 



d + 2 M 
2 W 



(33) 



According to Eqs. (|2"Tj) and ([3U]). the position at time t 
of a particle initially located at r = a is 



GAa d (l - Ba 2 )t. 



(34) 



The time at which the particle initially located at r 
a > reaches the origin r = is 



U(a) 



dGA(l - Ba 2 )' 



(35) 



In particular, the time at which the last particle, i.e. the 
one initially located at a = R, reaches the origin is 



2dGA dGM 



This is the final time since, at that time, all the particles 
have collapsed at the origin. The density profile is a 
Dirac peak p = M6(r) containing all the mass and there 
is no further evolution. As may appear surprising at first 
sight, the time at which the first particle, i.e. the one 
initially located at a — + , reaches the origin is finite: 



tcoll — 



2t]R d 



dGA d(d + 2)GM 



(37) 



This is the time at which the density profile becomes 
singular at the origin, since an infinitesimal number of 
particles has reached r = 0. This corresponds to a finite 



Combining the foregoing equations, the conservation of 
mass (HHJ) and the equation of motion (|2~T|) can be written 



J^a d (l 



M(r,t) = —^—a d (l-Ba 2 ), 



— (tcoll ~t)a d + B- — a 



d+2 



tr 



tr 



(39) 



(40) 



These equations, parameterized by a, completely deter- 
mine the evolution of the mass profile for all times. How- 
ever, in order to exhibit the self-similar structure of the 
solution close to t co ii (see Sec. IIVBI) . it is convenient to 
define 



1 a a 



tcou {tcoii - ty/ 2 



and 



in which case, Eqs. (|39| and (|40j) can be rewritten 



C = Bt 



(d+2)/d 
oil ' 



(41) 



(42) 



M(r,t) = j-(tcon t) d ' 2 y (l C t -^^y 2 / d ] , (43) 

"Lr \ t co ll J 



y + C 



tcoll 



(d+2)/d 



{tcoll -t)("+2)/2d' 



(44) 



(45) 



where x is the appropriate scaling variable and y is just 
a dummy variable. Using Eqs. ([Tg]), (|43]) and (02]), the 
corresponding density profile is 



(36) p(r,t) 



1 



SdG tcoii ~ t 



C 



2 + d 



1 — (tcoll — t) 

tcoii a 



,2/d 



The central density behaves like 

£ 1 



P (0,t) 



SdG tcoii — t 



(46) 



(47) 



and it diverges at t = t co ii- Finally, according to Eq. 
(|4"0"|) . the position at time t of the last particle, i.e. the 
one located at a — R at time t = 0, is 



r max (t) d = -^-(tcoii - t)R d + B-^R d + 2 . 

*coll *coll 



(48) 
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We can check that p(r max ,t) — 0. Let us now consider 
particular limits of these equations. 

For t — > 0, we obtain after some calculations 



B. Self-similar solution for t — > t co ii and r — > 



p(r,t)~p(r,0) + 



^col 



1 B 2(d+l)(d + 2) r2 



For t -> t co u and r -> 0, Eqs. (|43 |) -(|45 ]1 reduce to the 
form 



d 2 



_ B2 (d + 2)(d + 4)_ 4 



rf 2 



M(r,t) = ^(i coH -t) d / 2 2/ ) 



(49) 



The corresponding mass profile is 



M(r,t) ~ M(r,0) + 



_J ^ 

Gtcoll tcoll 



r d 2(d+l) d+2 

d d 2 



+S 2 



,d + 2 



d 2 



This is valid for r < r max (t) with 

2 i 



i 



d(d + 2) i coH 
For t — > 0, the central density behaves like 



p(0,t) 



SdGt 



coll 



1 



t 



t 



coll 



(50) 



(51) 



(52) 



Taking t = t co ii in Eqs. (|!3" ]) -([i5" |) . and using x d 
(jy{ d + 2 )/ d j we obtain the exact profiles 



M(r, t coU ) = 



1 



i - 



dG C<V( d + 2 ) 

(jd/(d+2 



r d 2 /(d+2) 



2d/{d+2) 



tcoll 



(53) 



1 



1 



SdGd + 2 C d /( d + 2 ) r 2d /( d + 2 ) 



d+2C d/(d+2) 

d £ co zi 



.2d/(d+2) 



(54) 



According to Eqs. (|48f and (|42[) . the position at time 
t co u of the last particle is 



C x l d 

r max — l d +2)/cP 
tcoll 



R (d+2)/d^ 



(55) 



We can check that M(r max ,t co u) = M. More precisely, 
using Eqs. (|5"3")) and (f5~4")) . we find that 



M(r,t coU )~M 



d :i 



2(d + 2) 



d 2 M 



1 



, (56) 



(57) 



for r r max . The behavior of M(r,t co u) and p(r,t co ii) 
for r — > is given in the following section. 



a; d = y + Cy 



{d+2)/d 



(tcoll -t)W>/™' 



(58) 



(59) 



(60) 



showing that the evolution becomes self-similar as we ap- 
proach the collapse time. The corresponding density pro- 
file is 



P(r,t) 



1 



1 



S d Gtcoii-tl+ d ± 2 Cy 2 / d ' 
The central density increases like 

p(0,t) 



(61) 



(62) 



and diverges at the collapse time. In parallel, the typical 
core radius 



r Q (t)~ (tcou -ty d+2) / 2d . 



(63) 



appearing in Eq. (|60|) decreases and tends to zero at t = 
tcoii- Comparing Eq. (|6"3"|) with Eq. (|S2"1) . we find that the 
core radius is related to the central density po(t) = p(0, t) 
by the relation 



Po^ ~ 1, 
involving the scaling exponent 

2d 



(64) 



(65) 



d+2 

Note also that the mass within the core radius is 

M(r (t), t) ~ p (t)r d (t) ~ (tcou - t) d/2 , (66) 

and it tends to zero as t — > t co ii- Therefore, there is no 

Dirac peak in the pre-collapse regime. 

For t = tcoii and r — > 0, Eq. (|59)l reduce to x d = 
Cy (d+2)/d and 

we obtain the power-law profiles 



M(r,t coU ) = 



1 



dG C d /^ 



1 



r d 2 /(d+2) 



1 



SdGd + 2 C rf /( d + 2 ) r 2d /( d + 2 ) ' 



(67) 



(68) 



p(r,t co u) 



They can also be derived from Eqs. (|53|) and ([54]) for 
r — > 0. This shows that the system develops a finite time 
singularity: the density profile diverges at the origin and 
scales like r - 2d /( d + 2 ) , By contrast, the mass M(r,t co ii) 
tends to zero at the origin. This means that only an 
infinitesimal number of particles has reached the origin 
at t = tc u, leading to infinite central density but zero 
central mass. 
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C. The half central density radius 

Using Eq. (|46|). the radius r*(t) corresponding to half 
the central density, i.e. such that p(r*(t),t) = p(0,t)/2, 
is 



r,(f) 



d 1\ 1/2 (tcoll - t)( d +^ 2d 



d + 2C 



x 2- 



t 



(d+2)/2d 
1/d 



d + 2t 



coll 



For t = 0, we obtain 



r*(0) 



R_ 

71' 



(69) 



(70) 



a result that can be directly obtained from Eq. ([31]) . For 
t — > 0, we find that 



2/ = 



-1 + Jl + AC^-x 2 



2Ct 



X 



trail — t 



(76) 



(77) 



For t — > t co u and r — > 0, the self-similar solution can be 
written 



M(r,t) = -L(t coll -t)y, 



P(r,t) 



J 1 1 

2nG t co u -t2Cy + V 



-1 + Vl + 4Cx 2 
2C : 



tcoll — t 



(78) 
(79) 

(80) 
(81) 



r.(t) 



71 



d 2 +4rf + 8 t 
4d(d + 2) i co/; 



(71) 



For t = tcoii, Eq. (p9^) gives = 0. In the self-similar 
regime t — > i co ;; . we find that 



.(*) = 



/ d 1\ 1/2 


'2(d+iy 


\d+2CJ 


d + 2 



1 1 / d 



(tcoll -tY d+ ^ 2d . 



(72) 



In the self-similar regime, the half central density radius 
r*(t) behaves like the typical core radius ro(i), see Eq. 
(|63|) . The mass within the sphere of radius r* evolves 
with time like 



M(r*(t),t) ~ PQ (t)r4t) d 
2 — 2 * 

d+2 ton u _i\d/2 
- (d+2)/2 ^ coU L > ■ 



2- 



t 

tcoll 



(73) 



In the self-similar regime, we find that M(r*(t),t) tends 
to zero like - (t co ii - i) d/2 - M(r (t),t). 



V. THE POST-COLLAPSE REGIME 

The previous solution shows that the system forms a fi- 
nite time singularity: the central density becomes infinite 
in a finite time t = t co ii- However, the mass contained 
within a sphere of radius e tends to zero as e — > 0. This 
is due to the fact that only an infinitesimal number of 
particles has reached r = at t = t co ii leading to infinite 
central density but zero central mass. This cannot be 
the final equilibrium state of the system since statistical 
mechanics predicts that the equilibrium state of a self- 
gravitating gas in the canonical ensemble is a Dirac peak 
containing all the particles (this is the configuration that 
makes the free energy tend to — oo due to the divergence 
of energy) [U, |45| . This is true at any temperature. For 
a self-gravitating Brownian gas a T = 0, we have previ- 
ously established that the Dirac peak p(r) = M5(r) forms 
at a time t en d given by Eq. (|36[) . Therefore, the evolu- 
tion continues in the post collapse regime t co ii < t < t en d 
where a Dirac peak grows and accretes progressively all 
the surrounding particles. The exact description of this 
process for T = is the object of the present section. 



D. More explicit solutions for d = 2 



In d = 2, Eq. |44|) can be solved explicitly to obtain ten 
y(x,t). Therefore, the exact solution can be written 



M(r,t) = ^(t coll -t)y 



tcoll 



(74) 



P(r,t) 



i 



2ttG tcoii — t 



, 2C , 

1 - - — (tcoii - t)y 

tcoll 



&+1 



(75) 



A. The exact solution 



For t > tcoii i the equation of motion (|21[) can be writ 



r d = a d - -GM(a, t co ii)(t - t coU ), (82) 

where M(a,t co ii) is the mass profile at time t = t co ii- 
Using Eq. (|53[) . we obtain 

1 



d d 

r = a — 



(Jd/(d+2) 



d 2 /(d+2) 



(jd/(d+2) 



2d/(d+2) 



t 



col! 



(t — tcoll)- 



(83) 
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At t = t c 



■At, the mass contained within the sphere of 



radius a* at t 
we get 



t co n has reached r = 0. From Eq. 



2d/(d+2) 



At 1 

C d/(d+2) i + AL 



(84) 



This leads to the formation of a Dirac peak of mass 
M D (t) = M(a*,t col i). Substituting Eq. flM]) in Eq. (|53)) . 
we obtain 



M D (t) 



dG C d l 2 



(d+2)/2 

(t~t c 



tr.nll 
t 



\d/2 



(85) 



We can check that Mn(t en d) = M so that all the mass 



has been absorbed in the Dirac peak at t = t e 
precisely, we find that 



M D (t) ~ M 



1 - 



2d 



1 - 



More 



(86) 



for t — > t enc i- It is relevant to write the conservation of 
mass (flU]) in the form 



M tot (r,t) = M D (t)+M(r,t) = M(a,t coU ), 



(87) 



where Mo{t) is the mass contained in the Dirac peak and 
M(r,t) is the mass exterior to the Dirac. Similarly, we 
write the density profile as 

p tot (r, t) = M D {t)5{v) + p(r, t). (88) 

The regular part of the profile p(r, t) is the residual den- 
sity defined as the density after the central peak has been 
subtracted. Using Eq. (|55|) . the evolution of the mass 
profile in the post-collapse regime is given by 



M D (t) + M(r 7 t) = 
x 1- 



1 



dG C d /<> d + 2 ) 

(jd/(d+2) 



7 d 2 /(d+2) 



d d 

r — a — 



1 



Qd/(d+2) 



2d/(d+2) 



J 2 /(d+2) 



(89) 



1 - 



(jd/(d+2) 



2d/(d+2) 



cull 



(t — t co n). 



(90) 



These equations completely determine the evolution of 
the mass profile in the post-collapse regime. In order to 
exhibit the self-similar structure of the solution close to 
t co u (see Sec. IVB[) . it is convenient to define 



f 



1 d 2 /(d+2) 



y = 



C d/(d+2) {t _ tcou) d/2> 

in which case, Eqs. (jg§)) and ([9"0)l can be rewritten 
M D (t) + M(r,t) - J-(t-t coll ) d/2 



dG 



xy[l-C- 



_tcM y 2/d 
tcoll 



(92) 



C—y 

tcoll 



(d+2)/d 



(i-t co n) (d+2)/2<r 



(93) 



(94) 



Note the similarities and the differences with the pre- 
collapse solution (|4"3")) - (|45[) . If we subtract the contribu- 
tion of the Dirac peak, using Eq. (1551) . we obtain 

M(r,t) = ^(t~t coll ) d / 2 



[1 _ C ttcoU y 2 /d 
tcoll 



, , \ (d+2)/2-| 

1 / tcoll \ 



Qd/2 



1 



x " = C—,/'+ 2 >' d - y, 

tcoll 



(i-t co /0 (d+2)/2d ' 



(95) 



(96) 



(97) 



Using Eqs. (TT51) . and (|M|) . the corresponding density 
profile is 



P(r,t) 



1 



SdG t — t co u 



l-^-{t-tc,u) 2 -^y 2 ' d 

tcoll « 

1 



d t co ii " 



(98) 



According to Eqs. dM]), & and (|H7]), the central resid- 
ual density decreases with time like 



p(0,t) 



i 



i 



SdG t — t co u V t 



tend 



- 1 



(99) 



Let us now consider particular limits of these equations. 

For t = tcoii, we recover the results (|53p - ([57| . For 
t = t end , using Eq. (|48]), we can check that r max (t en d) = 
0. More precisely, for t — > t en d, we obtain after some 
calculations 



P(r,t) 



dSdG tcoii 



te 



t 



tend 



G d/2 



d + 2\ d/2 r d 



(d+2)/2 



t 



'ill 



(100) 



The corresponding mass profile is 



(91) M(r,t) 



d 2 Gt 



coll 



t /d+2\ d/2 r d C d ' 2 



'-■ end 



2 A d + 2 )/ 2 

L coll 



(101) 



This is valid for r < r max (t) with 



Lend 



1/<1 



R. 



(102) 
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The total mass outside the Dirac is M(t) = M(r max ,t). 
Using Eqs. (jTUTj) and (fTD2]l . we obtain 



M(t) ~ M- 



2d \ t 



^end 



(103) 



Comparing with Eq. (JSBJ) , we check that Md (i) + M (t) = 
M. as it should. On the other hand, for t —> tend, the 
central density (|99|) decreases like 

p(0,t)~-^J-(l- J-V (104) 

in agreement with Eq. (| 100[) . 



B. Self-similar solution for t —> t co u and r — > 



appearing in Eq. f|108[) increases. Comparing Eq. (jllip 
with Eq. (|110p . wc find that the core radius is related to 
the central density po(t) = p(0,t) by the relation 



involving the scaling exponent 

2<7 



d+2' 



(112) 



(113) 



the same as in the pre-collapse regime. The mass within 
the core radius is 

M(r (t),t) ~p Q {t)r d {t) ~ (t-t coU ) d '\ (114) 



and it increases with time. Note that the mass (|105[) 
contained in the Dirac peak has a similar scaling. 



As the Dirac peak grows in the post-collapse regime, 
the central residual density decreases. Soon after t = 
t co u, the residual density profile p{r,t) follows a post- 
collapse self-similar evolution, reverse to the pre-collapse 
self-similar evolution, during which the central residual 
density p(0,t) decreases while the core radius ro[t) in- 
creases. This is the object of the present subsection. 

For t — > t co u , the mass of the Dirac grows like 



M D {t) 



dG C d / 2 



(t - t coll ) d/2 . 



(105) 



Furthermore, for t -> t coU and r -> 0, Eqs. ([55 ]) -([9T p 
reduce to the form 



M(r,t) = J^(t-t coll ) d / 2 



1 



£ = Cy 



(d+2)/d 



(t - t co «) (d+2)/2ti ' 



(106) 



(107) 



(108) 



showing that the evolution is self-similar just after the 
collapse time. The corresponding density profile is 



p(r,t) 



1 



S d Gt-t con ^Cy 2 /d-i 

The central residual density behaves like 

d£ 1 



p(0,t) 



2SdG t — t co u 



(109) 



(110) 



Therefore, the central residual density decreases as more 
and more mass is absorbed in the Dirac peak. In parallel, 
the typical core radius 



r (t)-(t-t coli r +2 ^ 2d , 



(111) 



C. The half central density radius 

In the self-similar regime, the core radius ro(t) in- 
creases with time indicating that the residual density 
profile expands as the central density decreases. How- 
ever, since the system size r max {t) decreases to zero for 
t = tend, the expansion of the density profile cannot be 
valid for all times. In fact, the core radius Tq{£) makes 
sense only during the self-similar regime. To study the 
evolution of the residual density profile in a more gen- 
eral setting, we introduce the half central density radius 
r,(t). 

Using Eq. (|9"5)l , the radius r*(t) corresponding to half 
the central density, i.e. such that p(r*(t),t) = p(0,t)/2, 



.(*) 



coll 



C d + 2 



1/2 



'1 , 2+d f„ou \V 2 
,2 ' 4 t ) 



4 ^coll T" 2 
1/d 



(d+2)/2d 



l l tcM -dh) ^- t ^ {d+2)/2d - ( 115 ) 

For t = t co u, this formula gives r* = 0. In the self-similar 
regime t — > t co u , we have 

/ s 1/2 / vl/rf 



C l/2 

For t — > t en d, we obtain 
rt{t) 



(116) 



1 



C d / 2 \d+2j It 



d/2 



tend — t (d+2)/2 



''end 



a . (H7) 



a result that can also be directly obtained from Eq. 
(|100[) . In the self-similar regime, we find that r*(t) oc 
(t — t co a)( d+2 >> 2d behaves like r (t) so that the density 
profile spreads. Then, ?%,(£) decreases and behaves like 
r* oc (tend — t) 1 '* for i — > i erl( j so that the density profile 
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shrinks. The mass within the sphere of radius r* evolves 
with time like 



M(r.(t),t) ~ Po{t)n{t) d 

(+ +\2 (1 I 2+dt coU \ d / 2 

V-end t>) \2 ' 4 t I 



t (d-2 



^ icoii i 



t s(d+2)/2 



(t -t coll ) d/2 . (118) 



In the self-similar regime, we find that M(r*(t),t) in- 
creases like ~ (t - t coU ) d/2 ~ M(r (*),*). Then, 
M{r*(t),t) decreases and behaves like M(r*(t),t) ~ 
{tend - t) 2 for t -> i 



D. More explicit solutions for d = 2 

In <i = 2, Eq. can be solved explicitly to obtain 

y(x,t). Therefore, the exact solution can be written 



/ ^ — teolt \ 1 / t co u 

tcoll ) C \ t 



2-i 



p(r,t) 



J 1_ 

2ttG t - t, 



coll 



2C , 

1 - - — \t - t co ti)y 

tcoll 

1 



(119) 



(120) 



VI. ILLUSTRATION OF THE RESULTS 

We shall now illustrate the previous analytical results 
by plotting some relevant quantities. For simplicity, we 
consider a two-dimensional system (d = 2) [65J and 
choose a system of units such that M = R = t co a = 1. 
With these conventions, we have A = 2, B = 1/2 and 
C = 1/2. On the other hand, £/G = 4 and S 2 = 2tt. 

The initial density profile (t = 0) is the parabole 



ir 



(127) 



The finite time singularity occurs at t = t co u = 1 and 
the Dirac peak containing the whole mass is formed at 
t = tend = 2. The size of the system, i.e. the radius at 
which the density vanishes, decreases like 



1/2 



(128) 



The exact density profile for < t < t co u = 1 (pre- 
collapse) is 



P(r,t) 



1 2 1 



1 - 1 7T VT+ 2tx 2 ' 



- — - { y/l + 2tx 2 - 1 



1 - t 



(129) 



(130) 



y = 



l + Jl + ACj^x 2 

~2Ct ' 

tcoll 



(121) 



For short times t — > 0, it is approximately given by 



p ( r ,t)~p(r,0) + ~ ^l-3r 2 + -r 4 )t, 



(131) 



t - tr 



For t — > t co u and r — > 0, the self-similar solution can be 
written 



M(r,t) = ^(t- tcou) (y-^ 



P(r,t) 



2nG t - t cM 2Cy-V 



1 + Vl+ACx 2 
2C : 



t-t c 



(123) 
(124) 

(125) 
(126) 



for r < TV 



1 — £/4 (see Fig. HJ. For < i < t, 



coll 



(122) the central density increases like (see Fig. [5]): 



, , 2 1 
p(0,t) = - : 



7T 1 - t 



(132) 



For t — > t co ;; = 1 and r — > 0, the density distribution 
takes the self-similar form 



Pss{r,t) = 



1 



i -r v l — t 



with the invariant profile 



i 



7T v/i + 2a; 2 



(133) 



(134) 



The central density corresponding to the self-similar so- 
lution ([TBTtj) is 



Pss(0,t) 



2 1 

ttI -t' 1 



(135) 
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and it exactly coincides with Eq. (|132|) . The singular 
profile corresponding to the self-similar solution (|133[) at 



t = t r 



1 is 



V2 



p ss (r,t coU ) = — . (136) 
irr 

The exact density profile at t = t co u = 1 is 

\/2 (1 



p{r, t co ii) = — ( - - V2 
7r \r 



(137) 



for r < r max = 1/ \/2. The exact density profile is plot- 
ted in Fig. [3] using self-similar variables. This repre- 
sentation illustrate the fact that the solution becomes 
self-similar for t — > t co u = 1; indeed the curves tends to 
the invariant profile f|134|) . The exact density profile is 
also plotted in Fig. 2] in logarithmic variables. This rep- 
resentation illustrates the fact that the density displays 
a finite time singularity at t — > t c0 u = 1 and that, ap- 
proaching the singularity, the profile becomes self-similar. 
For sufficiently small t co u — t, the tail of the profile (for 
ro(i) <C r <C r max (t)) is well-approximated by the pure 
power law (|136[) but if we want to describe the distribu- 
tion up to r max (t) ~ l/v^, we must use the truncated 
power law (|137j) . 




FIG. 1: Evolution of the density profile p(r, i) for small times. 
We have taken t = 0,0.05,0.1,0.2,0.3. The full lines corre- 
spond to the exact profile (|129l) and the dashed lines to the 
approximation (|131[) valid for i< 1. 



For t co ii = 1 < t < t en d = 2 (post-collapse), the mass 
contained in the Dirac peak grows like (see Fig. [5]): 



M D (t) = ^(t-l) 
The exact residual density profile is 



(138) 



p{r,t) 



1 



1 



t - 1 7T VI + 2tx 2 ' 



1 - 



t - 1 



1 + \/l + 2tx 2 



(139) 




FIG. 2: Evolution of the central density in the pre and post 
collapse regimes. The dashed line corresponds to the self- 
similar solution valid for t — > t co u — 1. 




FIG. 3: Evolution of the density profile p(r,t) us- 
ing self-similar variables. We have taken t = 
0, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 0.95, 0.99, 0.999. The 
full lines correspond to the exact profile (|129ll - (ll30[) . They 
tend to the invariant profile (| 134[) represented by a dashed 
line (hardly visible). 



(140) 



t-l 

For t — > t en d — 2, it is approximately given by 
2 (1-t 



P(r,t) 



(141) 



for r < Tr, 



- [(2-i)/2] 1 /2 ( S ee Fig. EJ). For t coU = 1 < 
t < t en d = 2, the central density decreases like (see Fig. 



, \ 2 1 /2 
p(0,t) = -— — --1 

7T t — 1 \E 



(142) 



For t — > t CO H = 1 + and r — > 0, the residual density dis- 
tribution takes the self-similar form 



Psa(r, t) 



t-l V t — 1 



(143) 
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FIG. 4: Evolution of the density profile p(r,t) in log-log plot. 
We have taken t coU -t = 1CT 6 , 1CT 5 , 10 -4 , 10 -3 , 10 -2 , 10 _1 . 
We have plotted the exact profile (|129|l - (|130|l and the self- 
similar profile (|133|1 - (|134[) (dashed lines). At t = t co u, the 
exact profiles tend to the density (|137p and the self-similar 
profiles to the density (|136[) . For sufficiently small t co u — 
t, the self-similar profiles and the exact profiles are in good 
agreement except for r ~ r max (t). 




0.08 - 



0.06 



°- 0.04 - 



0.02 




FIG. 6: Evolution of the density for t — > t end = 2. The 
full line corresponds to the exact profile and the dashed 
line to the approximation (|141|l . We have represented t — 
1.8, 1.9, 1.95, 1.97 (top to bottom). 




FIG. 5: Evolution of the mass contained in the Dirac peak. 
The dashed line corresponds to the self-similar solution. 



FIG. 7: Evolution of the density profile us- 
ing self-similar variables. We have taken t = 
1.001, 1.01, 1.05, 1.1, 1.2, 1.3, 1.4, 1.5, 1.6, 1.7, 1.8, 1.9, 2. The 
full lines correspond to the exact profile (|139[) - (|140[) . They 
tend to the invariant profile (|144|l represented as a dashed 
line. 



with 



/(*) = 



1 



7T Vl + 2X 2 



(144) 



The central residual density corresponding to the self- 
similar solution (|143[) is 



/ s 2 1 
Pss (0,t) = -- 



nt-1 



(145) 



The singular profile corresponding to the self-similar so- 
lution (|143j) at t = t co u — 1 is given by Eq. ()136|) and 
the exact density profile at t = t co u = 1 is given by Eq. 
(|137p . The exact density profile is plotted in Fig. [7] using 
self-similar variables. This representation illustrates the 
fact that the solution is self-similar for t — > t co u = 1 + ; 
indeed the curves tends to the invariant profile (|144p . 



the same as in the pre-collapse regime. Note, however, 
that the exact central density (|142[) is different from the 
self-similar central density (|145[) [they coincide only for 
t — > t co ii = 1 + ], contrary to the pre-collapse regime. This 
explains why the normalized central density (t— l)p(0, t) 
is not a fixed point in Fig. [7] contrary to Fig. [31 On 
the other hand, this figure shows the disappearance of 
the residual density profile at t = t enc [ = 2 when all the 
mass has been absorbed in the Dirac peak at r = 0. 
The exact density profile is also plotted in Fig. [5] in 
logarithmic variables. This representation illustrates the 
singularity at t = t co a and the fact that the solution is 
self-similar for t —> t co u = 1 + . As in the pre-collapse 
regime, for sufficiently small t — t co u, the tail the profile 
(for r (t) r <C r max (t)) is well-approximated by the 
pure power law (|136[) but if we want to describe the distri- 
bution up to r max {t) ~ l/y/2, we must use the truncated 
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1 




^ss^w Self-similar 




Exact I : 







10~ 6 1(T 4 10 2 10 ( 



r 

FIG. 8: Evolution of the density profile p(r,t) in log-log plot. 
We have taken t - t coU = 10~ 6 , 10 -5 , 10 -4 , 1CT 3 , 10 -2 , 10 _1 . 
We have plotted the exact profile (|139|l - (|140|l and the self- 
similar profile (|143|H144[) (dashed lines). At t = toll, the 
exact profile tend to the density (|137p and the self-similar 
profile to the density (I136|l . For sufficiently small t — t co u, the 
self-similar profiles and the exact profiles are in good agree- 
ment except for r ~ r max (t). We start to see a difference in 
the central density for t — t co ll = 0.1. This difference was im- 
perceptible in the pre-collapse regime for the equivalent time 
tcoii — t = 0.9. This is due to the difference in the laws of 
evolution of the central density (|132|l . (|135p [identical] and 
(114211 . (114511 [different] in the pre and post collapse regimes. 
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FIG. 9: Evolution of the radii r max (size of the system), ro 
(scaling variable) and r* (half-central density radius) with 
time. 



power law (|137|) . 

We represent in Fig. [9] the evolution of the half central 
density radius. For < t < t co u = 1. it is given by 

\-t ( t\ 1/2 
rS) = y t (2-{) , d46) 

and for t coU = 1 < t < t end = 2, by 

1 /2 

r,(t) = iQ + ij {2-tf'\t-l). (147) 



For t — > 0, we have 

T2 ('-!)• (148) 

and for t — > t en( i = 2, 

r*(t)^ l -(2-t) 1 ' 2 . (149) 
For t — > t co u = 1 j we find that 

r*(t)c±(-) \t-l\. (150) 

In the pre-collapse regime, the half central density radius 
decreases as the profile becomes more and more concen- 
trated. In the post-collapse regime, just after the col- 
lapse time t co u, the half central density radius increases 
as the profile expands, then reaches a maximum and fi- 
nally decreases to zero as all the particles are absorbed 
in the Dirac peak. In the self-similar regime, r*(t) has 
the same scaling as the core radius ro(t). 




t 

FIG. 10: Evolution of the moment of inertia / = (r 2 ) with 
time. 

Finally, we represent in Fig. [TU] the evolution of the 
moment of inertia 

rr mam (t) 

I = (r 2 ) = / p(r,t)r 2 2nrdr. (151) 

Jo 

It can be calculated from Eqs. (|129p and (|139|) after 
lengthy computations or, more directly, from the virial 
theorem (|B1[) of Appendix iBl For < t < t co u = 1, we 
find that 

I(t) =~{-\, (152) 
and for t co u = 1 < t < t en d = 2, we find that 

I(t) = ~(2-tf(3t-2). (153) 

Note that I(t) = (r 2 )(t) is a monotonically decreasing 
function of time. 
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VII. CONCLUSION 

In this paper, we have provided an exact analytical so- 
lution describing the gravitational collapse of a gas of self- 
gravitating Brownian particles in the overdamped limit 
at T = 0. Starting from a parabolic density profile, the 
system first develops a finite time singularity at t = t co u 
in the pre-collapse regime, followed by the formation and 
growth of a Dirac peak in the post-collapse regime. In- 
terestingly, our solution describes all the phases of the 
dynamics. This extends the purely self-similar solution 
obtained previously (36l - |38j that is valid only close to 
t co u. Unfortunately, our method that exploits the deter- 
ministic behavior of the system, is not valid anymore for 
T > 0. In that case, we must resort to other methods [36l — 
39| that are much more complicated. However, the phe- 
nomenology of the collapse at T ^ (finite time singu- 
larity, growth of a Dirac peak, self-similar solutions, pre 
and post collapse...) remains the same (but, of course, 
the scaling exponents are different). Therefore, for illus- 
tration of the general process of gravitational collapse of 
self-gravitating Brownian particles, it is useful to have a 
simple analytical solution such as the one described here. 

It is interesting to compare the isothermal collapse 
[H, EH of self-gravitating Brownian particles in contact 
with a heat bath to the collapse of isolated stellar sys- 
tems experiencing a gravothermal catastrophe (27l . [28| . 
These two types of systems presents analogies and differ- 
ences. For both systems, the pre-collapse is sclf-similar 
and generates a finite time singularity where the central 
density is infinite (see l36l . |37| for self-gravitating Brow- 
nian particles and [46l - l49| for stellar systems). However, 
for stellar systems, the post-collapse regime leads to a bi- 
nary star surrounded by a hot halo [46[ (statistical equi- 
librium state in the microcanonical ensemble [25l 12810 
instead of a Dirac peak containing all the particles [38j 
(statistical equilibrium state in the canonical ensemble 
44, HH). The binary star can release sufficient energy 
461 ] to stop the collapse and even drive a re-expansion 
of the system [5(| • Then, a series of gravothermal oscil- 
lations should follow [5l|. More references and discus- 
sions about the statistical mechanics of self-gravitating 
systems in microcanonical and canonical ensembles are 
given in the reviews [III [26| . 

For certain systems, such as those discussed in the In- 
troduction, it is important to take into account both long- 
range and short-range interactions. As we have seen, due 
to the attractive long-range interaction, self-gravitating 
Brownian particles, bacterial populations and colloids 
driven by attractive capillary interactions can collapse. 
In that case, the central part of the system becomes 
very dense. In the absence of short-range interactions, 
the collapse generically leads to the formation of Dirac 
peaks. Of course, in practice, Dirac peaks are unphys- 
ical and the density profile is regularized by small-scale 
constraints. These small-scale constraints can be due to 
finite size effects (the particles always have a finite size 
and cannot interpenetrate), stcric hindrance, short-range 



interactions and, ultimately, quantum mechanics (Pauli 
exclusion principle). These interactions come into play 
when the system becomes dense enough. Their effect is 
to provide a nonlinear pressure that can halt the col- 
lapse and lead to a well-defined equilibrium state. An 
example of this regularization is provided by a gas of 
self-gravitating fcrmions in which gravitational collapse 
is balanced by the pressure force arising from the Pauli 
exclusion principle I2ll [26l |52| . Another example is pro- 
vided by chemotaxis [12|| where Dirac peaks are replaced 
by smooth aggregates. In these examples, the Smolu- 
chowski equation at T = studied in this paper must 
be superseded by more general equations of the form (fT]) 
where the pressure p{p) prevents complete collapse of the 
system. A microscopic justification of these equations is 
given in 0, 0, [l3[ and some explicit examples of short- 
range regularization are worked out in [l2l . |52| . 

Finally, we would like to point out that the numerous 
analogies between self-gravitating Brownian particles, 
chemotaxis, colloids and some nanosystems described in 
the Introduction lead to the fascinating possibility of re- 
producing or mimicking a gravitational dynamics and a 
gravitational collapse in the laboratory. Experimental re- 
alizations of these processes will certainly be developed 
in the future and represent an interesting challenge. 



Appendix A: Nonlinear mean field Fokker-Planck 
equations 

The mean field drift-diffusion equation (fT]) is a subclass 
of the more general equation 

^ = V • {Dh(p)Vp + x.9(p)V$), (Al) 

coupled to Eq. ((SJ or to a potential $ = p * u (where 
u(r, r') is a binary potential of interaction and * denotes 
the convolution product) introduced in Here, both 

the diffusion coefficient Dh(p) and the mobility X9(p)/p 
can depend on the density. The fact that the diffusion 
and the mobility depend on the density in complex sys- 
tems is not surprising. For example, in a dense fluid, we 
understand that the motion of a given particle is ham- 
pered by interactions with its neighbors so that its mo- 
bility is reduced and its diffusivity modified (with re- 
spect to a dilute medium). For long-range interactions, 
these equations generalize the standard Debye-Hiickcl, 
Smoluchowski-Poisson and Kcller-Segel models. For 
short-range interactions, we can make the gradient ex- 
pansion $ ~ — ap — §Ap with a = —Sd J^°° u(q)q d ^ 1 dq 

and b = — ^Sd J^°° u(q)q d+1 dq [53|, and we obtain a 
generalization of the Cahn-Hilliard equations. Equation 
(|Aip was introduced in 0, Q as a nonlinear mean field 
Fokker-Planck equation associated with a generalized 
thermodynamical formalism. Later, Holm & Putkaradze 
[3^ considered a very related model 



dp 



V-(DVp + fi(p)pV<t>), 



(A2) 
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where $ = p * u and ~p = p * H, where u(\r — r'|) is a 
binary potential and H(\r — r'|) a local filter function. It 
differs from Eq. (|A1| only in the replacement of p by p 
in the diffusion and mobility terms. 



Appendix B: Virial theorem in d = 2 



In Appendix H of Ref. [54] . one of the authors has 
obtained the proper form of the virial theorem associated 
with the 2D Smoluchowski-Poisson system which is valid 
both in the pre and post collapse regimes. Its general 
expression is given by 



1 dl 

4** = " m 



k B T GM D (t) GM(t) 



(Bl) 



where I = J pr 2 dr = M(r 2 ) is the moment of inertia, 
Moit) the mass contained in the Dirac peak and M(t) = 
M — Mo{t) the residual mass contained in the regular 
density profile. In the absence of Dirac peak, Eq. (|B1[) 
reduces to the standard expression (55| : 



1 dl 

involving the critical temperature 



To), 



GMm 



This equation is closed and can be integrated into 



I(t) 



(T-T c )t + I a 



(B2) 



(B3) 



(B4) 



This is valid for any time when T > T c and in the pre- 
collapse regime when T < T c . In the post-collapse regime 
when T < T c , we must come back to the general expres- 
sion (|B1[) . However, this equation is not closed since it 
involves the mass Mb (t) contained in the Dirac peak. In- 
troducing the mean square displacement (r 2 ) = I(t)/M, 
Eq. (|B4I) shows that, for T > T c , the motion of a particle 
is diffusive with an effective diffusion coefficient j55j: 

^IH 1 -^)' (B5) 

that is independent on the initial condition. 
At T = 0, Eq. dBl|) reduces to the form 

tf t =-G[M 2 -M D (tf}. 

In the pre-collapse regime, t < t co u , we obtain 

GM 2 

I(t) = — t + I . 

In the post-collapse regime, t co u < t < t en d, we have 
established in Sec. [V] that the mass of the Dirac peak 
increases like (in d = 2): 



Substituting this expression in Eq. (|B6[) and integrating 
Eq. (|B6[) between t and t en d at which I(t en d) = 0, we 
find that 



1 ft 



I(t) 

MR 2 12 V t 



t 



coll 



3t 



coll 



(B9) 



We have checked, after lengthy calculations, that Eqs. 
(|B7[) and ()F39|) can also be obtained directly from the ex- 
act analytical density profiles (|75j) and (|120|) . The evo- 
lution of the moment of inertia (proportional to the dis- 
persion of the particles) is plotted in Fig. [TU] 

There is no closed expression of the virial theorem 
when d ^ 2. Therefore, the moment of inertia must be 
calculated directly from the exact analytical expressions 
of the density profile given in the main part of the pa- 
per. We shall give only particular values of the moment 
of inertia in d dimensions 

r r mam (t) 

I(t)= p(r,t)r 2 S d r d - l dr. (BIO) 

Jo 



For t ->• 0, using Eqs. (gSj) and ([BTjl. we find that 

I(i)~i- -_ — ^ MR 2 — . (Bll) 



For t = t co u, using Eqs. (|55]l and (|55]l. we get 



I (t coll) 



d :i 



2/d 



1 



-MR 2 . (B12) 



d + 2 \d + 2J d 2 + 2d + A 
For t — > t end , using Eqs. (|100[) and (|102j) . we obtain 



2{d? + 3d + 2)" 



end 



2(rf+l)/d 



(B13) 



Using Eqs. (|102[) and (|103p . we check that this last ex- 
pression can be written 



J(t) 



rf 2 + 3d + 2 



(t) 2 . (B14) 



For d = 2, the expressions (jBlip . (|B12j) and (|B13|) are 
consistent with the general results (|B7I) and (|B9[) valid 
for all times. For the uniform sphere in d dimension, 
using Eqs. P5 ]) -([2^ ]| . we obtain 



(B6) m 



- — -Mi? 2 1 - — — j-t 
d + 2 V 



2/d 



d + 2 



MR{tf 



(B15) 



(B7) In d = 2, this expression reduces to Eq. (|B7[) as it should 



Appendix C: Derivation of the generalized mean 
field Smoluchowski equation 



M D {t)=AM{ t ?f) (-^-1 

t I \ tcoll 



(B8) 



In this Appendix, following [l3j, we provide a new 
derivation of the generalized mean field Smoluchowski 
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equation JT]). In previous works @, H, [l(J, this equa- 
tion was derived from a notion of generalized thermo- 
dynamics. In that case, the nonlinear pressure was due 
to a bias in the transition probabilities leading to non- 
Boltzmannian distributions @. Here, we show that the 
same equation can be derived from the Dynamic Density 
Functional Theory (DDFT) used in the theory of simple 
liquids [56]. In that case, the nonlinear pressure is due to 
the correlations induced by the short-range interactions. 
Although physically distinct, these processes lead to the 
same type of macroscopic equations. A more detailed dis- 
cussion, as well as the derivation of more general kinetic 
equations, is given in fl3( | . 

We consider the overdamped dynamics of N Brownian 
particles in interaction governed by the coupled stochas- 
tic equations Q: 



* =-,V (! 7( ri ,.. 



,r N ) + V2DHi(t), 



(CI) 



where U(ri, r N ) = m 2 J2i<j u (\ r i ~ r j\) i s the poten- 
tial of interaction and Ri(i) is a Gaussian white noise 
such that (Ri(t)) = Oand (Rf{t)R^(t')) = S i:j S a pS(t-t'). 
Here, i = 1, N label the particles and a = 1, d the 
coordinates of space. The diffusion coefficient D is re- 
lated to the mobility p = 1/(£to) and the temperature 
T by the Einstein relation D = pk B T The time 
evolution of the iV-body distribution Pjv(*"i, ...,rjv,t) is 
governed by the TV-body Fokkcr-Planck equation 



N 



dt 



N „ 

Y — 

^dr z 



k B T dP, 



N 



m 



dr. 



1 9 
-P N —U(n, 



,rjv) 



(C2) 

This particular Fokkcr-Planck equation is called the JV- 
body Smoluchowski equation. Its steady state is the 
Gibbs canonical distribution Pjv = -^e~^ u . 

It is easy to derive from Eq. (|C2[) the equivalent of the 
BBGKY hierarchy for the reduced distribution functions 
Q. Introducing the local density p(r,t) = NmPi(r,t) 
and the two-body distribution function ^(r, r',i) = 
N(N—l)m 2 P 2 (r, r',t), the first equation of the BBGKY- 
like hierarchy is the exact Smoluchowski equation 



dp 
dt 



k R T. 



X7p+ p 2 (r,r',t)Vu{\r-r'\)dr' 



(C3) 



where we have used the fact that the particles are iden- 
tical. This equation is not closed since it involves 
the two-body correlation function p2(r,r',t). We must 
therefore introduce some approximations to evaluate this 
term. We shall assume that the potential of interac- 
tion u = ulr + use is the sum of a long-range poten- 
tial ulr. and a short-range potential usr- For systems 
with long-range interactions, it is known that the mean 
field approximation is exact in a proper thermodynamic 
limit N — > +oo Therefore, concerning the long- 

range potential u_l_r, we shall make the approximation 



p 2 (r,r',t) = p(r,t)p(r',t) leading to 

J p2(T,T',t)WULR(\v-r'\)dv' = p(T,t)V$(T,t), (C4) 



with 

$(r, t) = J p(r', t)u LR (\r ~ r'|) dr' . (C5) 

To evaluate the integral corresponding to the short-range 
interactions, we shall use an approximation that has be- 
come standard in the dynamic density functional theory 
(DDFT) of fluids [1^ and take 

SF 

p 2 (r, r', t)Vu SR (\r - r'|) dr' ~ p(r, t)V°-^[p(r, t)], 

(C6) 

where F ea; [p] is the excess free energy [66| calculated at 
equilibrium. This relation is exact at equilibrium (58j 
and the approximation consists in extending it out-of- 
equilibrium with the actual density p(r, t) calculated at 
each time. This closure is equivalent to assuming that 
the two-body dynamic correlations are the same as those 
in an equilibrium fluid with the same one body density 
profile. Although it is not possible to ascertain the valid- 
ity of this approximation in the general case, it has been 
observed for the systems considered that this approxi- 
mation gives remarkable agreement with direct Brownian 
TV-body simulations. With the approximations (|C4|) and 
~2>J), Eq. dOS]) becomes 



dp 



ksT „ „5F ex 
-2-Vp + pV— — 

m op 



(C7) 



which is closed. The total free energy is 

F[p] = I [ p$ dr + k B T [ In ^dr + F ex [p] . (C8) 
2 J J m m 



The steady state of Eq. (|C7|) minimizes (|C8|) at fixed 
mass and is given by SF—pSM = where p is a Lagrange 
multiplier. This yields SF/Sp = p i.e. 



p(r) = Ae 



(C9) 



In a fluid, the local pressure is of the form p = p(p, T). 
Since the temperature T is fixed in the case of Brow- 
nian particles (canonical description), the pressure is 
barotropic and we shall simply write p = p{p). In prin- 
ciple, the excess free energy -F^p] can depend on the 
gradients of the density. This is particularly important 
for a fluid close to an interface j58|. Here, we shall as- 
sume that the density varies on a distance that is large 
with respect to the range of intermolecular forces. This 
is the case if the density distribution is mainly due to 
long-range interactions, as we shall assume in the follow- 
ing. With this assumption, the free energy is of the form 



F\P\ = \ J P$>dr + J p f^kp-dp! dr. 



(CIO) 
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The excess free energy is therefore 



F ex [p] 



Pi 



We note the relation 



(Cll) 



k T SF 
Vp(p) = —Vp + pV—^=V Pld + V Pex , (C12) 
m dp 

where Pid(r) = p(r)fcsT/m is the ideal gas law and p ex 
the excess pressure due to short-range interactions [59|. 
This relation can be used to determine the equation of 
state p(p) corresponding to the excess free energy F ex [p] 
and vice versa. For an ideal fluid (F ex = 0), we recover 
the perfect gas lawp(r) = p(r)kBT/m. For a free energy 
of the form ([UlOj) . using identity ([CT2"]) . Eq. ([C7]) can be 
rewritten 



pV$). 



(C13) 



The equilibrium state is given by the condition of hydro- 
static equilibrium 



generalized mean field Smoluchowski equation (|C15|) iff 
it is a critical point of F at fixed mass. Furthermore, 
it can be shown that a steady state is linearly dynami- 
cally stable with respect to the generalized Smoluchowski 
equatio n (|C15[) iff it is a (local) minimum of F at fixed 
mass 0, 3| . This is consistent with the condition of ther- 
modynamical equilibrium. If F is bounded from below 
[67| . we know from Lyapunov's direct method that the 
system will converge towards a (local) minimum of F at 
fixed mass M for t — > +oo. If several (local) minima 
exist (metastable states) , the choice of the selected equi- 
librium will depend on a complicated notion of basin of 
attraction. 

Remark: note that Eq. (|CJ15|) can be justified in a phe- 
nomenological manner from the linear thermodynamics 
of Onsager if we interpret it as a continuity equation 
dtp + V • J = with a current J = -{l/£)V5F/5p pro- 
portional to the gradient of a potential p(v) = SF/Sp 
that is uniform at equilibrium (see Eq. (|C9I0 . This is 
precisely the way in which this equation was introduced 
in the physics of liquids (58| and, more generally, in 0, 0] . 



Vp + pV$ = 0. 



Equation (|C13|) is a generalized mean field Smoluchowski 
equation including a gcnerically nonlinear barotropic 
pressure p(p) due to short-range interactions and a mean 
field potential $ = Ulr* p due to long-range interac- 
tions. As we have indicated in the Introduction, this 
equation arises in several physical problems such as self- 
gravitating Brownian particles [3, chemotaxis (2{|, col- 
loids with capillary interactions [30(, etc. Combining re- 
sults issued from the physics of systems with long-range 
interactions and from the dynamic density functional 
theory of fluids [56j , this equation has been justified from 
a microscopic model [HI . 

The generalized mean field Smoluchowski equation 
(|C7[) can be written in terms of the free energy functional 
£30) as 



dp 
dt 



1 &F 

i sp 



(C14) Appendix D: The formation of a peripheric Dirac 

peak 

In Sec. we have mentioned that the evolution of 
the system is qualitatively different whether M(a, 0)/a d 
is an increasing or a decreasing function of a. In the main 
part of the paper, we have considered the physical situ- 
ation where M(a,0)/a d decreases. We have shown that 
such initial conditions lead to the growth of a Dirac peak 
at the origin. In this Appendix, we illustrate on a spe- 
cific example the situation where M(a, 0)/a d increases. 
We show that such initial conditions lead to the forma- 
tion of a peripheric Dirac peak (a d-dimensional annu- 
lus) that progressively converges towards the center of 
the domain by absorbing the interior particles. To sim- 
plify the formulae, we choose a system of units such that 
M = £ = G = 1. We also assume that the initial size of 
the system is R = 1. 
(C15) An analytical solution can be obtained by considering 

an initial condition of the form 



This equation monotonically decreases the free energy 
functional (|C8[) which plays therefore the role of a Lya- 
punov functional. Indeed, a straightforward calculation 
leads to the ff-theorem appropriate to the canonical en- 
semble 



F = 



SFdp , 

Jp-m dr = 



1 ^ SF \ 



SF 

V— ) dr < 0. (C16) 
op. 



For a steady state, F = 0, the last term in parenthesis 
must vanish so that SF/dp is uniform. This leads to Eq. 
(IC9[) . Therefore, a density p(r) is a steady state of the 



M(a,0) = a 



2d 



(Dl) 



for a < 1 and M(a, 0) = 1 for a > 1. The corresponding 
density profile, given by Eq. (fT%]). is 



p(a,0)=™a d , 



(D2) 



for a < 1 and p(a, 0) = for a > 1. The density increases 
with the radius for a < 1 and vanishes discontinuously 
at a = 1. According to Eq. ([2Tjl. the position at time t 
of the particle initially located at a is 



da 2d t 



(D3) 



18 



According to Eqs. (|19[) and (|D1I) . the mass profile at 
time t is 



M (r, t) = a 



2</ 



(D4) 



where a is related to r and t through Eq. (|D3|) . The 
size of the system at time t, i.e. the position of the last 
particle initially located at a — 1, is 



r max (t) = {l~dt) l t d . 



(D5) 



From Eq. (|D5[) . all the particles have collapsed at r = 
at the final time 



1 

tend ~j ■ 



(D6) 



Solving for a in Eq. (|D3| , which is a second degree equa- 
tion in a d , and substituting the resulting expression in 
Eq. (|D4| . we obtain 



(D7) 



According to Eq. (|T5|) . the corresponding density profile 



is 



1 



1 



- 1 



(D8) 



It is represented in Fig. [TT] at different times. Since 
t < Tmax (t) , we can check that the term under the square 
root is always positive. It vanishes for r = r max at the 
collapse time 



1 

tcoii - 



(D9) 



At that time, the density becomes infinite at r = 
r m ax(t C oii) = (l/2) 1 ' d . On the other hand, we note that 
the central density p(0, t) = at any time. The density 
profile is self-similar since 



p(r,t) = - t f(rt^ d ), 



with the invariant profile 



fix) 



S d Wl - 4:dx d 



(D10) 



(Dll) 



This function increases and diverges when x — > (l/4(i) 1 / d . 
We must distinguish two cases: 

(i) For t < tcoii-, a simple calculation shows that 



p(r m ax(t),t) 



2d 1 

5^1 - 2dt' 



M(r max (t),t) = 1. 



(D12) 



(D13) 




FIG. 11: Evolution of the density profile p(r,t) in d = 2. The 
time t goes from to t en d = 0.5 and is represented every 0.05 
time steps. The vertical line corresponds to the peripheric 
Dirac peak, appearing at t = t co n — 1/4 and containing all 
the mass at t — t en d = 1/2. 




FIG. 12: Evolution of the mass contained in the peripheric 
Dirac peak in d = 2. 



In that case, the profile (1D7|) - (|D8[) contains all the mass 
and there is no Dirac peak. The density at the periph- 
(t),t) increases monotonically with time and 
diverges when t = t co ii ■ The density profile at t — t co u is 



p(r, t C ou) 



2d 



1 



S d V Vl - 2r d 



- 1 



(D14) 



for r < r max (t coU ) = (l/2) 1 / d and p(r,t) = otherwise, 
(ii) For t > t co u, a simple calculation shows that 



(D15) 



(D16) 



In that case, the profile (|D7[) - (|D8|) docs not contain all the 
mass. A Dirac peak forms at r = r max (t) and captures 
the peripheric mass as the system shrinks. The mass 



p{r max {t),t) = 



(t),t) 



2 1-dt 
~S~{i2dt-l' 



(dty 
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contained in the Dirac peak is Mu(t) = 1— M(r max {t), t) 
yielding 



M D (t) 



2dt-l 



(D17) 



It obviously satisfies Mrj{t co u) = and Moitend) = 1 
(see Fig. [T2")) . The density at the periphery is infinite at 
t = t co ii and monotonically decreases for t > t co u as the 
Dirac peak grows by absorbing the neighboring particles. 
The total (normalized) density profile can be written 



Ptot(r,t) = p(r,t)+M D (t) 



5{r 



m) 



S d r 



d-1 
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where the first term p(r, t) is the regular density profile 
(|D8j) and the second term corresponds to the Dirac peak. 
For t —t t co u and r — > r max {t)~ , the regular density 
profile (|D8[) has the self-similar form 



p{r,t) 



1 
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We can understand the onset of the Dirac peak for- 
mation in a more qualitative manner. According to Eq. 
(|D3|) . the positions of the particles initially located at a\ 
and <22 > a\ coincide at the time 



1 

t = -- 



1 



and take the value 



n = r 2 



of 



(D20) 
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The shortest time at which this coincidence occurs cor- 
responds to ai ~ ci2 ~ 1 yielding t = t co u = l/(2d) 
and n = ?'2 = (l/2) 1 / d . At that time, the density be- 
comes infinite since there is a finite mass in the inter- 
val r 2 — r\ — > 0. Then, a peripheric Dirac peak with 
mass Mi){t) forms at r max (t) and grows by capturing 
the neighboring particles as the system collapses. At 
t = tend-, the Dirac annulus has reached the origin and 
all the particles are contained in a central Dirac peak at 
r = 0. 

Remark: According to Eq. (|D3[) , we would naively 
conclude that the particle initially at a reaches the origin 
r = at time t*(a) = l/(da d ). This time decreases 
with a so that the particles that are initially close to the 
origin seem to take more time to collapse at r = than 
far away particles. In particular, for a < 1, we find that 
t*(a) > t enc i which leads to an apparent paradox. In 
fact, the above estimate is not correct. Indeed, a particle 
initially located at a is captured by the peripheric Dirac 
peak at a time t c (a) such that r(t c (a)) = r max (t c (a)). 
This leads to t c (a) = (1 - a d )/[d(l - a 2d )]. For a -> 1, 
we find that t c (a) t co u = l/(2d) which corresponds to 
the onset of the peripheric Dirac peak. We also note that 
t c (a) < t <t {a) for a < 1, so that the particles are captured 



by the peripheric Dirac peak before they can reach the 
origin according to their pure motion (|D3|) . Therefore, 
they are carried towards the origin by the Dirac peak 
whose equation of motion is given by Eq. (|D5I) . This 
solves the apparent paradox. 



Appendix E: Repulsive interaction 

In the main part of the paper, we have considered the 
situation where the interaction between particles is at- 
tractive like gravity (G > 0). In this Appendix, we 
briefly consider the case of a repulsive interaction like in 
a Coulombian plasma (G < 0). In chemotaxis, when the 
chemical substance secreted by the biological entities acts 
as a pheromone, the interaction is attractive (chemoat- 
traction) but when it acts as a poison, the interaction 
is repulsive (chemorepulsion) . The case of a repulsive 
interaction can therefore have interesting physical appli- 
cations. 

In the deterministic case (T = 0), the general exact 
solution of the problem is given by Eqs. (jT9")l and (f2"Tj) by 
simply making the substitution G — > — G. Accordingly, 
the dynamical evolution of a uniform sphere is given by 
Eqs. ([2"5 ]) -(f2"9" )l with G replaced by -G. The sphere re- 
mains uniform during the evolution but it now expands 
instead of collapsing. The mass profile, the radius and 
the density evolve like 



M(r, i) 
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r 



R d 1 



dGM 



t 



dGM 



l/d 



P(t) 



P(0) 



dGM 4. ' 

eR 3 1 
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(E3) 



For t —> +oo, we obtain the scalings M(r, t) ~ £r d / (dGt), 
R(t) - (dGMt/tl) 1 ^ and p{t) - £/{S d Gt). 

For the initial parabolic profile (|30|) - (|33|) . the exact 
time dependent solution can be written 



M(r,t) = Aa d (l - Bo?), 



GAa d {\ - Bo 2 ) 



r 

x =—d- 
The corresponding density profile is 

F 1 - 



p(r,t) 
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(E7) 



dGAt 
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The size of the system, i.e. the position of the last particle 
initially located at a = R at time t = 0, is 



,(t) = R[ 1 + 



dGM 
£R d 



l/d 



(E8) 



We note that this expression coincides with Eq. (|E2[) for 
all times. The central density (r = a = x = 0) decreases 
like 



p(0,t) 
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S d Gti 
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dGAt 



For t +oo, the system behaves like a homogeneous 
sphere since Eqs. (|E4|) - (|E9|) become equivalent to Eqs. 
(|E1|) - (|E3[) in this limit. Considering the general equations 
of the problem, we note that this property is valid for an 
arbitrary initial profile. For t — > 0, we find that 



p(r,t) ~ p(r,0) 
1 



d 2 {d + 2) 2 GM 2 
4^ Jr™ 
2(d+l)r 2 d + 4r 4 



d R 2 d+2R 4 
The corresponding mass profile is 



t. 



(E10) 



M(r,t) ~ M(r,0) 



d(d + 2) 2 GM 2 d 



£,R 2d 
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d + 2 i? 2 d + 2i? 4 _ 
These expressions are valid for r < r max {t) with 
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0, the central density evolves like 
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In d = 2, Eq. (]E5[) can be solved explicitly to obtain 
a(x,t). We get 
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2B 



2GAt 
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Substituting this expression in Eqs. (|E4|) and (|E7[) . we 
obtain an explicit solution for M(r, i) and p(r, t). To sim- 
plify the notations, we introduce dimensionless variables 
such that M = R = 1, A = 2, B = 1/2, £/ G = 4 and 
& = 2ir like in Sec. ED Equations (|B7} and (|E14j) can 
be combined to give 



1 



a/(1 + i) 2 - 2r 2 t 



(E15) 




FIG. 13: Evolution of the density profile p(r, t) in d = 2 for 
t = 0, 0.5, 1, 2, 3, 5, 10. The initial density profile is parabolic 
but it asymptotically tends towards a step function (uniform 
disk). 



for r < r max (t) with 



:(*)= 1 + 



1/2 



(E16) 



For r > r max (t) , we have p{r, t) = as usual. The density 
profile is plotted in Fig. [13] for different times. For t = 0, 
the profile is parabolic and for t —> +oo, it tends to a 
Heaviside function (uniform disk). The central density 
decreases like 



P(0,t) 



2 1 



nl + t 
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In the repulsive case, making the substitution G — > — G 
in Eqs. (|B2|) - ()B3|) . the virial theorem in d = 2 dimensions 
reads 



ldl / GMm 
-^= N ( kBT+ — 
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It can be integrated into 



J(t) = — U B T+— ^— )t + J(0). 



For T = 0, we obtain 



CM 2 

J(t) = — t + J(0). 



(E20) 



Introducing the mean square displacement (r 2 ) = 
/(t) /M, Eq. (|E19[) shows that the motion of a particle is 
diffusive with an effective diffusion coefficient 



D(T) 



k B T GM 



(E21) 



that is independent on the initial condition. For T = 0, 
we obtain Dq = GM/(4£). For the uniform sphere in d 
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dimension, using Eqs. (|E1|) - (|E3|) . we get 

d , / dGM \ 2/d d , v9 

(E22) 

In d = 2. this expression reduces to Eq. (|E20[) as it 
should. For a parabolic initial condition, introducing the 



dimcnsionlcss variables defined previously and using Eq. 
(IE15|) . we find after some calculations that 

lit) = (r 2 ) =\ + \- (E23) 
This expression coincides with Eq. (|E20[) . as it should. 
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We can obtain a larger class of models (see [j| and Ap- 
pendix^ by letting the friction coefficient £ depend on 
the density, £ = £(p), or more generally be a function of 
position and time, £ = £(r, t). 

To avoid misunderstandings, we stress that the 
Smoluchowski-Poisson system does not describe tradi- 
tional astrophysical systems. Indeed, astrophysical sys- 
tems are not in a strong friction limit but rather in a weak 
(or no) friction limit. For example, stars are described by 



22 



the Euler-Poisson system (hydrodynamics) galax- 
ies by the Vlasov-Poisson system (collisionless dynamics) 
[22l | and globular clusters by the Vlasov-Landau-Poisson 
system (collisional dynamics) (2^|. The Smoluchowski- 
Poisson system could describe the dynamics of dust par- 
ticles in the protoplanetary nebula, where particles ex- 
perience a friction with the gas and a stochastic force 
due to turbulence or other diffusive effects [24j. Unfor- 
tunately, when the dust particles are small, the friction 
force is strong but the gravitational interaction is weak. 
Alternatively, when the particles are large, the gravita- 
tional interaction is strong but the friction force is weak. 
Therefore, the Smoluchowski-Poisson system is only valid 
in an intermediate regime which may not be the most 
relevant. Nevertheless, in the context of the statistical 
mechanics of self-gravitating systems [25I. [2(1] . the model 
of self-gravitating Brownian particles is interesting at a 
conceptual level because it provides a dynamical model 
associated with the canonical ensemble [jO] while stellar 
systems are associated with the microcanonical ensemble 
[27I . [2^ | . It can be used therefore to illustrate the notion 
of ensemble inequivalence that is generic for systems with 
long-range interactions [f|. 
[63] The theoretical interest of these equations was stressed 
early in Q| before all their possible applications were re- 



alized. 

[64] The case of repulsive interactions G < is treated in 
Appendix [Ej 

[65] Note that the dimension d — 2 is particularly relevant in 
chemotax is [2911 and for the dynamics of colloids at fluid 
interface [30t | . 

[66] The excess free energy -Feajp] is a non-trivial functional 
determined by the short-range interactions. All the diffi- 
culty in the theory of fluids is to find some approximate 
forms of this functional. Once this functional is deter- 
mined, the density profile, as well as all the n-point corre- 
lation functions, can be obtained via functional differen- 
tiation. Inversely, the excess free energy is often obtained 
from the study of the correlation functions. The excess 
free energy F ex is known exactly only in a few particular 
cases, but very good approximations can be devised in 
more general cases [58LI59]. 

[67] This is not always the case. For example, the free 
energy associated with the Smoluchowski-Poisson sys- 
tem describing self-gravitating Brownian particles is not 
bounded from below [H, |37], 0, In that case, the 
system can experience an isothermal collapse. However, 
there also exists long-lived metastable states on which 
the system can settle [jjfl, 0, [gO] . 



